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Abstract 

A new three-dimensional model of the FEL is presented. 
A system of scaled, coupled Maxwell-Lorentz equations 
are derived in the paraxial limit. A minimal number of 
limiting assumptions are made and the equations are not 
averaged in the longitudinal direction of common radi- 
ation/electron beam propagation, allowing the effects of 
coherent spontaneous emission and non-localised electron 
propagation to be modelled. The equations are solved nu- 
merically using a parallel Fourier split-step method. 

INTRODUCTION 

A three-dimensional model of a helical wiggler Free 
Electron Laser (FEL) is presented that minimises the 
assumptions made. This model does not average the 
Maxwell-Lorentz equations describing the interaction be- 
tween the electrons, wiggler and radiation fields. Further- 
more, no relativistic approximations in the equations gov- 
erning electron motion are made and transverse motion of 
the electrons is self-consistently driven by both the wiggler 
and radiation fields. Current 3D codes perform averaging 
over a radiation wavelength of the radiation and others also 
over the the electron trajectories. Thus no study of the ef- 
fects of coherent spontaneous emission (CSE) or self am- 
plified coherent spontaneous emission (SACSE) is possi- 
ble. Furthermore, these codes cannot model large energy 
losses by the electrons as the electrons are confined locally 
to a wavelength. In this model no such averaging is done 
and so sub-period radiation evolution is included and elec- 
tron migration over distances greater than a wavelength, for 
example due to locally high radiation powers, may be de- 
scribed. 

3D MODEL 

Following the method of previous studies for the one- 
dimensional FEL LL 2J the physics of the FEL may be 
described in three dimensions by the coupled Maxwell- 
Lorentz equations 
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where j — 1 . . . N, the total number of electrons, and the 
transverse current density may be written 
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where m is the electron rest mass. The helical wiggler and 
radiation electric field are assumed to be 
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where e = (ic + iy)/^/2, ^0 is a slowly varying com- 
plex field envelope and B^, is the wiggler magnetic field 
strength of period A^, = 27r/fcuj. The magnetic compo- 
nent of the radiation field B(r, t) is required in the Lorentz 
equation and may be calculated from Maxwell's equations 
to give: 
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where it has been assumed that any radiation that is 
counter-propagating the electron beam may be neglected. 
The following scaUng notation is now introduced: 
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where (3^ — v^o /c is the mean electron velocity on entering 
the interaction region, 7^ is the resonant electron energy 
and p is the fundamental FEL parameter, defined as 
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Here ujp — ^ e'^Up/ (som) is the plasma frequency for the 
peak electron number density of the electron pulse Up, and 
Qw = cByj / (mckw) is the wiggler deflection parameter. 
The scaled electron density is defined as fip = Uplgl^, 
where Ig is the gain length and Ic the cooperation length. 
Note the scaling in the transverse plane is with respect to 
the physically relevant gain and cooperation lengths. This 
allows the equations above describing the FEL interaction 
to be written as 
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In deriving the equations the only approximations made are 
the neglect of space charge and the paraxial approximation. 
There are no restrictions on the electron energy allowing 
large changes to be modelled. 

In the scaled variables used here, the normalised beam 
emittance is given by 

En — \/lglcOxOp^ ■ (13) 

For an electron beam matched to a focussing channel the 
radius of the electron beam envelope is given by 
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The electrons experience a natural focussing force due to 
the magnetic field of the helical wiggler |4| resulting in a 
beta-function of /3 = This restoring force is 

approximated by the final linear term in (O and may be 
simply modified in magnitude to approximate other, e.g. 
FODO, focussing channels. 

The scaled Rayleigh range (in units of z) can be shown 
to be zr — (T?/2p, where CTp is the radius of the beam in 
scaled units of x,y. 

It can be shown that the above equations satisfy energy 
conservation, written in the form: 
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where 
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which may be used as a check in the numerical solution. 

In the relativistic Compton limit e, p ^ 1 the above 
equations can be shown to reduce to those of |1|| . 



NUMERICAL SOLUTION 

The field evolution of the scaled equations (l7li-(fT2]) are 
solved by applying a modified parallel Fourier split-step 
method Q in conjunction with a finite element method 
1^. This allows the effects of coherent spontaneous emis- 
sion (CSE), self-amplified spontaneous emission (SASE) 
and diffraction to be modelled numerically. The paral- 
lel code was developed using the Numerical Algorithms 
Group (NAG) parallel routines to control the parallel pro- 
cessing. The Fourier split-step method involves solving the 
wave equation in two steps. The first step considers the 
field diffracting in the transverse direction and propagating 
in the z direction freely without the electron source term: 
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This equation can be solved by taking the Fourier trans- 
forms in X, y and Z2 resulting in an ordinary differential 
equation in the Fourier transformed field A, with the fol- 
lowing analytic solution: 



A{z + Az) = A{z) exp 
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where the transverse Fourier transform variable pairs are 
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given by (x, fc^), (y, ky), {z2, fcjj and k\ 
inverse numerical Fourier transform is then applied giving 
solution A{x, y, Z2, z + Az) which is then used as the ini- 
tial field for the second part of the split-step method where 
the diffraction terms are neglected and the source term acts 
alone: 
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The finite element method is then used to solve the wave 
equations ( fTSl l and (l8Tl-(fT2]i along with a fourth order 
Runge-Kutta method for electron variables to give the fi- 
nal solution. 

Before applying the finite element method to the wave 
equation, the summation over the real electrons has to be 
changed to a summation over a group of macroparticles, 
each of which represents many real electrons. This reduces 
the computational memory load and operation: 
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where subscripts j and k indicate evaluation at the particle 
position, Np is the total number of macroparticles and Nk 
is the charge weight of the macroparticle in units of the 
electron charge. 

The Galerkin method of the finite elements is used. The 
field is described by a set of 8 -node hexahedral elements 



with linear basis functions Am{x,y,Z2) and nodal values 
am{z): 

A{x, y, z, Z2) = ^ a„(z)A„(5, y, Z2). (20) 

where m is the global node index of the 3D system. Note 
that Am(a;, y, Z2) = outwith the elements to which it be- 
longs. The wave equation ( fTsT i then reduces to a matrix 
equation for the nodal points: 

dz ~a^f^^ Ve (l + |pi|2)l/2 
I^m{xk,yk,z2k)e~ 

where K is the stiffness matrix constructed from the ele- 
mental equations f6l, % is a macroparticle weighting func- 
tion and Ve is the scaled volume of an element of the sys- 
tem, Ve — Aa;AyAz2. 

A SIMULATION IN THE ID LIMIT 

To test the 3D model a comparison with the ID results 
of m 121 was made. In these one-dimensional models the 
wave equation was integrated analytically over the com- 
mon electron/radiation transverse area before numerical so- 
lutions were obtained. 

Here, one element is used in the transverse plane approx- 
imating the radiation field to a plane wave and allowing 
integration over the transverse area and comparison with 
the ID results. The electron beam radius was set greater 
than the electron orbital radius ensuring the electron beam 
does not move significantly in the transverse plane of the 
element. A gaussian distribution for the electrons in the 
transverse plane was used with a range of six times the 
standard deviation. The results of Ul Fig. 2 showing self- 
amplified coherent spontaneous emission (SACSE) for a 
top-hat pulse current are reproduced in Fig. [T] (Note that 
there is a scale reversal as 1 1 1 plots the power as a func- 
tion of zi while here it is plotted as a function of Z2-) 
The coherent spontaneous emission, including the effects 
of shot noise, from both square and gaussian shaped pulses 
of ||2i were also reproduced. In Fig. |2] the top hat case 
of ||2| Fig. 4 is plotted. Good agreement between the ID 
models of |[Tl|2l and the ID limit of the 3D model of this 
paper are observed. 

A SIMULATION IN 3D 

An example of the code operating in 3D is now pre- 
sented. The parameters used are not intended to model any 
specific system, but rather to demonstrate the functioning 
of the code and its associated post-processing and plotting 
routines. A relatively simple system was chosen with a 
short electron pulse so that coherent spontaneous effects 
can be observed and large computational effort is not re- 
quired. 
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Figure 1: The scaled radiation power \A\'^ as a function 
of Z2 demonstrating good agreement with the ID model of 
HI Fig. 2 
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Figure 2: The scaled radiation power |Ap as a function 
of Z2 demonstrating good agreement with the ID model 
ofl2JFig.4 



The energy spread and the emittance of the beam were 
set to zero, the emittance by setting ap^ = for all elec- 
trons. Focussing was therefore not required for this simu- 
lation. The code was run using 16-processors of a parallel 
computer with parameters shown in Table [T] 

Fig.[3]plots the scaled intensity \ A\'^ across a transverse 
slice at the head of the electron pulse, just as the radia- 
tion escapes to propagate into vacuum. In Fig.|4]the scaled 
power (|/lp integrated across the transverse plane) is plot- 
ted as a function of Z2 following a FEL interaction to z = 5. 
Both plots demonstrate the main features of FEL interac- 
tion with Coherent Spontaneous Emission and the energy 
conservation relation of ( fTSl l is confirmed to within 0.5%. 

Further testing of the code is being carried out and will 
be reported upon in a future publication. 
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Propagation distance z 


5.00 


Radiation parameters 




Initial seed field over pulse Aq 


0.01 


Seed distribution in x & y 


Gaussian 


Sigma in a^.y 


~ 0.228 


Rayleigh length zji 


1.74 


Distribution of seed in Z2 


Top-hat 


Length of seed in Z2 


1.5 



Table 1 : 3D parameters 




Figure 3; The scaled radiation intensity |Ap in the trans- 
verse (a;, y) plane at z = 5. The transverse shce is taken at 
the 'head' of the electron pulse at Z2 = 5. 



CONCLUSION 

A new parallel code has been developed which mod- 
els the FEL amplifier by solving a system of scaled equa- 
tions describing the FEL interaction in three spatial and the 
time dimension. The aim has been to introduce as few 
approximations into the model as possible, the main as- 
sumptions being the neglect of space charge and any back- 
ward propagating radiation fields. This allows the effects of 
coherent spontaneous emission, diffraction and full elec- 
tron transport throughout the region of integration to be 
modelled. Furthermore, the sub-wavelength discretisation 




Figure 4: The scaled radiation power at z = 5 as a function 
of Z2- The electron pulse is shown on the Z2 axis lying 
between 5 < Z2 < 6.5. In this frame the electron pulse 
propagates along Z2 with increasing z. 

of the model allows a significantly wider range of radia- 
tion frequencies to be modelled than is possible with most 
other codes that use a minimum discretisation interval of 
a radiation wavelength. A ID limit of the computational 
model was identified and simulation results in this limit 
show agreement with previous ID numerical and analytical 
models. A further example simulation has been presented 
which demonstrates the code operating successfully using 
an electron pulse in three dimensions. Further optimisa- 
tion of the code is ongoing. While the code is undoubtedly 
slower to operate than other 3D averaged codes, the ex- 
tended physics that it can model may be expected to yield 
interesting new phenomena in FEL physics. 
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